Frozen shuffle update 
for an asymmetric exclusion process 
^ ; with open boundary conditions 

1-^ , C. Appert-Rolland, J. Cividini, and H.J. Hilhorst 

0^ , 1 - University Paris-Sud; Laboratory of Theoretical Physics 

Batiment 210, F-91405 ORSAY Cedex, France. 
2 - CNRS; UMR 8627; LPT 

^ ; Batiment 210, F-91405 ORSAY Cedex, France. 

2< . E-mail; Cecile.Appert-Rolland@th.u-psud.fr, 

Julian . CividiniSth . u-psud . f r , Henk . HilhorstSth . u-psud . f r 



C/3 



Abstract. We introduce a new update algorithm for exclusion processes, more suit- 
able for the modeling of pedestrian traffic. Pedestrians are modeled as hard-core par- 
C^ ' tides hopping on a discrete lattice, and are updated in a fixed order, determined by a 

phase attached to each pedestrian. While the case of periodic boundary conditions was 

studied in a companion paper, we consider here the case of open boundary conditions. 

rj ' The full phase diagram is predicted analytically and exhibits a transition between a 

O . free flow phase and a jammed phase. The density profile is predicted in the frame of 

a domain wall theory, and compared to Monte Carlo simulations, in particular in the 
vicinity of the transition. 
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1. Introduction 

As part of the effort to describe pedestrian motion, tfie past years fiave seen the 
development of cellular automata based models, among which the so-called 
'floor field' model [Ill2]- These models represent pedestrians as particles that 
jump from site to site on a discrete lattice, with an exclusion principle that forbids 
two pedestrians to simultaneously occupy the same site. In general, the pedestrians' 
positions are subject to a parallel update procedure, i.e., the particles attempt to jump 
at the same instants of time. Although this type of update ensures a quite regular 
motion under free flow conditions, it creates conflicts when two or more pedestrians try 
to move at the same time to the same target site. 

Introducing an arbitrary numerical decision procedure to resolve these conflicts can 
be considered as a drawback and alternative update procedures have been looked for. 
In reference [1], for example, large pedestrian simulations are performed by means of 
a random shuffle update. In this sequential update scheme an updating order for the 
pedestrians is drawn at random at the beginning of each time step. 

The random shuffle update was first proposed and studied in [5l [6], [7j. It was 
applied there to the Totally Asymmetric Simple Exclusion Process (TASEP), which can 
be seen as a basic model for pedestrian traffic. More generally, exclusion processes are 
archetypes of out-of-equilibrium systems and have been used to model various transport 
phenomena, including road and intracellular traffic, or ant motion |8j. 

Here and in a companion paper [Qj we introduce a new variant of the shuffle 
update called the frozen shuffle update. Under this updating procedure each particle 
i is assigned, either at the beginning of the simulation or when it enters the system, 
a random phase Ti G [0,1). In each time step particles are updated in the order of 
increasing phases. The phases may be thought of as the phases in the walking cycle, 
where one pedestrian may be are slightly in advance with respect to another one. This 
advance then determines who passes first in case of conflict. 

We have recently studied the frozen shuffle update for the TASEP on a ring, i.e. 
on a closed system [9]. This system is appropriately described by its current J as a 
function of the particle density p. We found [9] that the TASEP on a ring undergoes a 
phase transition at a critical density p = |, where J{p) attains its maximum through a 
linear cusp. The critical point separates a 'free flow' phase from a 'jammed' phase. 

In the present paper we study an open-ended TASEP with frozen shuffle update. 
Whenever a new particle, say i, enters through the open boundary, its place in the 
updating order is specifled by a newly drawn Tj. It will then keep this same ti until it 
leaves the system. In section [2] we deflne the frozen shuffle update. We will consider 
only its deterministic version in which allowed jumps are carried out with probability 
one. In addition to the bulk rules of motion, which are the same as on the ring [9], 
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we specify in section [3] the rules by which particles enter and leave the system. This 
procedure will depend on two additional parameters, an entrance probability a and an 
exit probability /3. In section H] the phase diagram of the system in the a/3 plane will 
be determined analytically. We find a transition line a = f3 between, again, a free flow 
and a jammed phase. The bulk particle density p{a,/3) and the current J{a,(3) are 
determined analytically and appear in excellent agreement with simulations. In section 
[5] the finite size effects, in particular near the transition line, are described by a domain 
wall approach. 

2. Frozen shuffle update for the ID TASEP 

2.1. Frozen shuffle update 

We consider a system of hard-core particles {i.e. at most one per lattice site) on a 
one-dimensional lattice of L sites numbered from left to right by k = 1,2, ... ,L. The 
particles enter the system at A; = 1, make hops of a single lattice distance to the right, 
and leave the system at k = L. The hops are executed according to the following frozen 
shuffle update scheme. A particle i, when it enters the system, is assigned a random 
phase Ti G [0, 1), which it keeps until it leaves the system. The phases are drawn from a 
predefined distribution. In each time step all particles are updated {i.e. attempt to hop 
one step to the right) in the order of increasing phases. An attempted hop is successful 
if and only if the target site is empty. 

Throughout, the variable t will stand for continuous time; time steps will be 
indicated by integer values t = 0,1, . . . ,s, . . .. The particle configuration at time t = s 
is considered to be the result of the sth time step. 

2.2. Two stable bulk states 

In recent work [9\ we have identified two distinct stable bulk states sustained by this 
update scheme. In an infinite homogeneous system these may be discussed without 
reference to boundary conditions. They are the free flow state and the jammed state. 
Free flow state. In the free flow state all particles are able to move forward at each 
time step. Hence they all have a velocity v = 1, that is, one lattice distance per time 
step. The current J^ = pv in the free flow state is therefore in these units given by 

J^ = p. (1) 

The technical definition is that the system is in a free flow state if whenever two adjacent 
lattice sites are occupied by particles, the particle to the right has the lower phase. 

Jammed state. The second stable bulk state is the jammed state. This is a highly 
compacted state of a special kind. Technically, the system is in the jammed state if 
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Figure 1. Six hard-core particles on a one-dimensional lattice whose sites are 
represented by square boxes. The arrow indicates the hopping direction. Our 
convention will be to place a heavy (red) vertical bar to the left of each particle 
that has a larger phase than its follower (in the present example, particle 
3). The configuration shown is not a free flow state, since the particle pair 
(3, 4) violates the free flow condition; particle 4 is unable to move during the 
upcoming time step. The conflguration is not a jammed state either, since 
the holes between the particle pairs (1,2) or (4,5) violate the jammed state 
condition. 



whenever a particle has a smaller (larger) phase than its follower, then it is located on 
the site immediately to the right of its follower (is separated from its follower by at most 
a single hole). As a consequence, in this state, a fraction of the particles remains blocked 
at each time step; the average particle velocity therefore is less than unity. 

Stability. The conditions for the free flow state and the jammed state are illustrated 
by figure [H It is tacitly understood that these conditions are to be evaluated after 
completion of a full time step. There exist many particle configurations that are neither 
free flow nor jammed state^jl. For a system on a ring lattice the latter two were shown 
[H] to be attractors in phase space toward which arbitrary initial states evolve. They 
are stable in the sense that once the system has entered one of them, it will stay in it. 

Platoons. We now look for the counterpart in the jammed phase of the current- 
density relation ([T]). The flow in the jammed state is most easily discussed in terms of 
platoons. A platoon is a succession of particles occupying consecutive sites and having 
increasing phases; the first particle of a platoon having a lower phase than the particle 
that precedes it; and the last particle of a platoon having a higher phase than the 
particle that follows it. It was shown [9] that in the jammed phase each particle belongs 

:|: There also exists a small class of configurations that are, technically, both 'free flow' and 'jammed', 
namely the configurations in which each platoon is preceded by a hole; this fact is of little consequence 
for the discussion. 
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to a platoon. 

Frozen shuffle update leads to simple platoon dynamics. A platoon moves as a 
whole, i.e., if its first particle can move, then in the same time step the whole platoon 
will move one lattice distance to the right. Equivalently, one may say that the hole 
initially in front of the platoon has hopped across it to the left over a distance equal 
to the platoon length (as will be further illustrated in section [3^ . This allows us to 
make the following statement [9]. Let u be the average platoon length in some suitable 
spatially homogeneous statistical ensemble. Since the hole density is 1 — p, and since 
only platoons that have a hole in front of them can move, the current J'^ in a jammed 
phase is given by 

JJ = ^(1 - p) , (2) 

which is the relation that we looked for. Henceforth the superscripts J and F will 
indicate quantities referring to the jammed and the free flow state, respectively. 

As the final comment of this section we remark that expressions ([I]) and (I2l) refer 
to spatially homogeneous states of, in principle, infinite extent. Our interest below will 
be in a finite system with boundary conditions specified by two parameters a and /3. In 
the next section we define these boundary conditions and express the stationary state 
current J{a, (3) of this finite system in terms of the currents J^ and J"^ discussed above. 

3. Rules for entering and exiting 

3.1. Entering rule 

The rule governing the entering of particles into the system is motivated by an equivalent 
hard rod description [9] of the system in the free flow statq§|: Let nonoverlapping rods of 
unit length, as shown in figure [2], move in continuous time at constant speed v = 1 along 
the X axis. At any integer instant of time s, each rod i will cover exactly one integer 
spatial coordinate x = k. We will say that at time step s the particle i associated with 
the rod i occupies lattice site k in the corresponding discrete model. This particle i will 
then be on lattice site fc + 1 at time s + 1, and so on. The (noninteger) time s + Tj at 
which rod i uncovers position k and starts to cover position k + 1 defines the phase r, 
of particle i in the discrete model. 

If the rods are (at any instant of time) uniformly distributed in space, subject 
only to the no-overlap constraint, then the intervals Tj between them (see figure [2]) are 
independent variables drawn from an exponential distribution 

P{T) = ae""^, T > 0. (3) 

where a is a positive parameter. 

§ This equivalence breaks down when there is jamming, as it would lead to overlapping rods. 
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Figure 2. Rods move at constant speed v = \ from left to right along the x axis. 
The time interval between the tail of the ith and the head of the (i + l)th rod 
is called Tj . The Tj are drawn randomly from the exponential distribution ([3]). 
We have also represented the underlying discrete lattice, that allows to define 
the mapping between the continuous and the discrete model. If the snapshot is 
taken at an integer time s, then the phase Tj of particle i in the discrete model 
corresponds in the continuous model to the distance between the tail of rod i 
and the integer position k that rod i overlaps (since u = 1, time and space may 
be identified). 



The equivalence between the rod and the particle motion lead us to adopt the 
following rule for particles entering the lattice, in order to have the same time interval 
distribution between particles in the discrete system. A particle can be injected from 
the exterior onto the first lattice site (labeled 1) only if that site is empty. If site 1 gets 
emptied at some instant of time s + rj_i (due to particle i — \ moving forward to site 2), 
then it will be reoccupied by a new particle i at a random time tj = s + rj_i + Tj with 
Tj distributed according to ([3]). This new particle will therefore have a phase Tj given 
by the recurrence relation 



Ti = tj mod 1 



(4) 



= (ri_i + Ti) mod 1, 

to which we will return in section [3J 

Below it will sometimes be convenient to use the probability, to be called a, that 
the initial site, when empty at some arbitrary (real) instant of time t, is occupied at 
time t + 1. This quantity is independent of t and with ([3]) it follows that 

a = Prob{T < 1} 

= 1 - e-^ (5) 

We will refer to a and a as the 'entrance rate' and 'entrance probability', respectively. 



Frozen shuffle update for open TASEP 7 

3.2. Entering current 

Let p^ and J^ denote the spatial density and the current, respectively, of the incoming 
particles under free flow conditions. They then satisfy the free flow relation ([T]), that is, 
J^ = p^ . We now express these quantities in terms of the parameter a. Let T^ be the 
time needed for A^ particles i = 1,2, ... ,N to enter the system. Then T^ is the sum of 
the time intervals Tj and, for each particle, of a unit time interval during which it blocks 
the access of a new particle to site 1. In the limit A^ — ?■ oo the stochastic variable T™ 
becomes sharply peaked around its average, that is, 

A'' 
rnm 1 1 

hm ±if = hm - y (Ti + 1) = 7l+ 1 = - + L (6) 

Since J^ = limiv_j,oo N/T^^ , it follows that at the left hand end of the lattice we have, 
under free flow entrance conditions, 

l+a ^ l+a ^ ' 

The injection procedure seeks therefore to impose onto the system a free flow state with 
a density and current ([7]) characterized by the parameter a (or a). This free flow is 
interrupted only when a particle meets with an obstacle: that might be a barrier at the 
exit point or a jam caused by other particles anywhere in the interior of the system. In 
the special case that site 1 is occupied by a particle which is itself blocked, the entering 
rule of section Ism still applies, but in that case the free inflow stops. Before commenting 
further we will now examine what happens at the exit point. 

3.3. Exiting rule 

We choose to adopt the following rule for particles leaving the system. If the last lattice 
site, labeled L, is occupied by a particle i with phase Tj, then, at each instant of time 
s + Tj, particle i leaves the system with probability /? and stays on site L with probability 
1 — /3. For /3 = 1 the free flow is unimpeded at the exit. However, when /3 < 1, the free 
flow will at times be randomly impeded, which may create a slowly advancing waiting 
line for the particles near the exit. When the rear end of this waiting line propagates 
all the way backward to the entrance, the system is in its jammed state. Below we will 
investigate the conditions for this to happen. 

3.J^. Exiting current 

Let the particles near the exit point be in a jammed state configuration and let p^ and 
J'' denote their density and their current, respectively, under this condition. These two 
quantities then satisfy the jammed state relation ([2]), that is, J'^ = v{l — p'^). We will 
now show how they are determined by the parameters /3 and a. 
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Figure 3. The lattice near the exit site L, shown at five successive time steps (s 
is an arbitrary integer). Particles in the same platoon are labeled by the same 
letter; their phases, not indicated, increase towards the left. A heavy (red) 
bar indicates an end-of-platoon. The particles are in a jammed state (platoons 
are separated by at most a single hole). The particle motion is deterministic 
except for the exit from site L, which takes place only with probability /3. Time 
step t = s: Initial configuration. Time step t = s + 1: The particle heading 
platoon a exits and its follower hops to the exit site. A new hole, marked by an 
up-arrow '|", enters the system from the right. Besides, platoon h advances by 
a position exchange with the hole marked by a cross x. Platoon c is blocked. 
Time step t = s + 2: Particle a, subject to the random exit rule, stays on site L. 
Platoons h and c advance. The two marked holes exchange positions with the 
platoon to their left. Time step t = s + 3: The last particle of platoon a leaves 
the system and is replaced with a new hole on site L. Platoon b is blocked but 
c and d move ahead. Time step ,3 = ^-1-4: Platoons h and d move ahead, but c 
is blocked. 
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We reason as follows. A particle on site L has at each time step a probability f3 to 
leave the system. The probability q^ that it exit at the rth time step is q,. = (1 — /3)'"~^/3, 
so that for leaving the system it needs an average number of time steps equal to 
YlT=i''^1r = t^~^- III the jammed state the particle configuration consists of platoons. 
An example is shown in figure [3l where the last particle of each platoon is indicated by 
a heavy (red) vertical bar to its left. 

When the particle on site L exits, the whole platoon that it was heading advances 
one step, so that site L is occupied by the next particle in the exiting platoon. When 
the last particle of a platoon has left the system, site L is unoccupied and remains so 
during one time step, after which it is occupied by the first particle of the next platoon. 
Therefore, the average number of time steps needed for an ra-particle platoon to exit, 
counted from the moment that its first particle arrives on site L until the moment that 
the next platoon's first particle arrives there, is equal to n/3~^ + 1. If z/ is the average 
platoon length, then z//3~^ + 1 is the average time needed for u particles to exit the 
system. At the right hand end of the lattice we therefore have, under jammed flow exit 
conditions, 

Before combining our analyses of the entrance and the exit points we complete this 
section by finding an expression for the mean platoon length u, appearing in ([8]). It 
should be noticed that z/ is completely determined by the injection procedure. Actually 
we shall see that it can be expressed as a function of a. 



3.5. Expression for z/(a) 



The average platoon length v depends on the injection procedure defined in section [3?T| 
and in particular on the way the phases Tj are assigned. Let tq be some arbitrary initial 
phase. The phases Tj are then related to the time intervals Tj by 
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where we set 



Tj = Tj modi. 
The probability distribution P for the Tj may be calculated as 
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= J2 dTae-^^dif-T + m) 

oo 

= J2 ae-'^^'^+™)^(f + m), (11) 

m=—oo 

in which 6{x) is the Heaviside step function, here defined by 6{x) = for x < and 
6{x) = 1 for a; > 0. Substituting for 6{x) its definition in flTTl) restricts the sum, which 
can then be calculated as follows 



Pm = J2 



^g-a(T+m) 



m=0 



a 



^-al ^ < T < 1. (12) 



1-e- 
A useful intermediate step is to rewrite ([9]) as 

Ti = fi mod 1 (13) 

with 

i 

f, = To + J2Tr (14) 

i=i 
The fj are a sequence of random points on the positive f axis, not confined to [0,1), 
and there is at least one point in every interval between two integers k and k + 1. If 
fi and fj+i are in the same interval [k, k + 1), then the corresponding Tj and Tj+i verify 
Ti < Tj+i after the modulo in ( IT3l) . If fj lies in [k, k + 1) and fj+i in [k + l,k + 2), then 
recalling that Tj+i < 1, we have Tj > Ti+i. Thus each set of fj values within an interval 
[k,k + 1) will give one platoon in the jammed phase. The average platoon length u is 
therefore equal to the average number of fj values in a unit interval on the f axis, which 
in turn is the inverse of the average length of the modulo 1 time intervals Tj. From flT2|) 
we now obtain easily 
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which is the desired expression for z/ in terms of a, or equivalently, a. 

Upon combining (IT5l) with ([8]) we find that the exiting current J'^ may be expressed 

as 
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We remark that rewritings like the one above are carried out most efficiently by working 
with a and a as though they were unrelated variables. 

This completes the necessary preliminaries. Everything is in place now for us to 
determine the phase diagram of the system subject to the boundary conditions {a, /3). 

4. Phase diagram 

In the previous sections, we have identified and characterized two possible stationary 
states, free flow and jammed. It is clear that the free flow state cannot be sustained if 
the outflow predicted by ( 1T6|) is smaller than the inflow given by (ITj); in that case the 
jamming will propagate across the system from exit to entrance. Inversely, it is also 
clear that the jammed flow state requires the outflow, given by ( IT6l) . to be smaller than 
the inflow ([7]). Hence in the stationary state viewed as a function of a and /3 we expect 
the occurrence of a phase transition on the line 

J^{a) = j\a,f3), (17) 

where we have indicated explicitly the dependence of the currents on the two model 
parameters a and /3. It is interesting to notice that ( jTBj) may be rewritten as 

Then one immediately finds the critical line J^ = J'' in the a/3 plane is 

a = 13. (19) 

This simple result would have been rather unsurprising in the case of particle/hole 
symmetry. However, in the present case the particle/hole symmetry is broken by the 
update scheme, and we have not found a simpler derivation of flT9|) than the one given 
above. 

As a consequence of (fT9|) we have for the current J{a, (3) carried by the system in 
its stationary state 

f J^(a), a</3, 
J{a,P) = { (20) 

I Jj(a,/3), a>/3, 

and 

p{a,P) = { (21) 

I p-^(a,/3), a> P, 

with J^ and p^ given by (JTj) and J"^ and p"' = J'^ / (3 by ( IT6l) . Across the critical line the 
current is continuous and the particle density jumps by a factor (3, 

J! = Ji, pi = PpI (22) 
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where the subscript 'c' indicates evaluation on the critical line. This jump is a true 
discontinuity only in the limit L — )■ oo. Finite size rounding will be described in the 
next section. 

The closing remarks of this section concern the similarities and differences that 
appear here between the present frozen shuffle update scheme and other schemes applied 
to the TASEP. First, we note that the critical line ( IT9l) obtained here for the TASEP is 
identical to the one found with other update schemes, whether it be random sequential 
update [lOl [H], parallel update [12], or ordered sequential update [I3]. Secondly, 
however, our expressions for J^ and J^^ are specific to the update procedure adopted 
here. In particular, it is remarkable that in the present case the jammed phase current 
depends not only on /3 but also, through the mean platoon length z/(a;), on the entrance 
probability a. Indeed, the jammed phase dynamics is completely determined by the 
platoon structure, which in turn results from the injection procedure - whatever this 
procedure may be. 

Thirdly, the frozen shuffle update of this paper, unlike the other updates, does not 
lead to a 'maximal current' phase: here the current is always equal to either J^ or J'^, 
that is, always fixed by the boundary conditions. This feature, confirmed by simulations, 
is due to the deterministic character of the TASEP considered here: when a particle has 
the possibility to hop forward, it does so with probability 1. 

5. Domain wall picture 

5.1. Microscopic domain wall position 

li (3 < 1, a waiting line may be formed at the exit whose rear end is nothing but a 
domain wall separating a spatial region with free flow from one with jammed flow. For 
a < /3 the system is in the free flow phase; the incursions of the jammed region into 
the system, starting from the exit point, will not exceed some finite localization length 
and be short-lived. Hence the domain wall is localized near the exit point. For a > (3 
the system is in the jammed phase; the waiting line invades the whole system and the 
domain wall is localized near the point of entry. For a ~ /3 the domain wall position 
is subject to large fluctuations that we will study below; for a = /3 it is completely 
delocalized. 

The frozen shuffle update of the TASEP offers the advantage, compared to 
previously studied update schemes, that the wall position may be defined very 
accurately, namely on a microscopic scalq]]]. Since in the free flow state no particle 
is blocked at any time, it is natural to say that the domain wall position is determined 
by the leftmost particle ever to have been blocked during an attempted move. By 

II A similar microscopic definition could be given for a parallel update with deterministic dynamics. 
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Figure 4. Density profile p{k^ + i) in the coordinate system attached to the 
domain wall; I is the distance to the fluctuating wall position k^. The data 
were taken on the critical line with a = (3 = 0.4. The profile is very close to a 
step function that jumps from p'l to p^ at the vertical (red) line. Only in the 
jammed state does there appear some structure on the scale of a few lattice 
distances. 



convention we will take for the domain wall position k^ the half-integer coordinate 
of the link on the left of this leftmost blocked particle. When the leftmost blocked 
particle moves again, it carries the wall along; when another particle positioned to its 
left is prevented to move, the wall position is transferred to the link on the left of that 
particle. All particles to the left of the wall position are in the free flow state and all 
those to its right are in the jammed state, in the sense of the definitions of section 12.21 
There is no need to develop tracking strategies as was necessary [H], for example, for 
the TASEP with random sequential update. 

Let p{k) be the site dependent particle density. Figure H] shows our simulation 
results for the profile in the reference system of the wall, within an environment of 20 
lattice distances around k^. Data were taken only when the wall was more than 20 
lattice units away from both ends. By construction p{k^ + |) = 1. Apart from this 
exceptional site the density profile is almost a step function, separating fiat density 
profiles with densities p^ and p'^ , as shown by the vertical (red) line in figure HI Only 
on the high density side of the wall does there appear some structure, restricted to a 
microscopic distance of the order of a platoon length (equal to z/(0.4) = 2.18 in the case 
of figure H]) . 
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Figure 5. Space averaged particle density p versus entrance probability a for 
the open-ended chain with exit probability j3 = 0.4. The simulation data are 
for a sequence of increasing system sizes. They converge towards the lower 
(black) and upper (red) solid curves, which are the theoretical densities p^ {a) 
and p^{a, (3) of equations (|31|) and (j32|) . The dotted curves are their metastable 
extensions. At the critical value a = ac = 0.4 the density jumps from its free 
flow value p'l to its jammed flow value p^. The simulation data for finite system 
sizes L, however, show strong finite size rounding. 



5.2. Domain wall theory 

The strong fluctuations of the wall position near criticality give rise to large finite size 
effects visible on the physical observables. One of these is the spatially averaged density, 
still denoted p, in a finite system. Simulation data for this quantity are shown in figure 
[5] and are seen to converge only slowly to the discontinuous theoretical curve. Below, a 
description of the system dynamics in terms of the wall motion will allow us to account 
for such finite size effects. 

Let Qk{t) be the probability for the wall to be located at k^ = k + ^, i.e. on link 
[k, k+1), where we may take k = 0,1, . . . ,L. The "links" (0, 1) and {L, L+1) are virtual; 
Qoit) is the probability that the leftmost blocked particle occupies position k = 1 and 
Qiit) is the probability that none of the particles present in the system at time t has 
ever been blocked. We will follow here the mesoscopic approach to domain wall motion 
that was earlier applied with success to the TASEP with other update schemes [151 [IB] 
and that has been able to predict stationary and nonstationary [T71 [HI [19] features of 
those TASEP versions. According to this approach we hypothesize that 
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(i) the probability Qk{t) obeys the master equation for an asymmetric simple 
random walk, 

^Qfc(t) = D+ [Qk-i{t) - Qk{t)] + D^ [Qk+i{t) - Qk{t)] , (23) 

for k = 0,1,2, ... ,L with the reflecting boundary conditions Z^+Q-i = D-Qo and 
D^Ql+i = D+Ql, and 

(ii) the coefficients D± in this equation are related to the currents and the densities 
in the two phases by 

^+ = rr^' D- = ^^- (24) 

P'' — p^ p-' — p^ 

For the present case of frozen shuffie update the validity of equations (!23|) and ( 124|) has 

not been demonstrated on the basis of first principle^Sl- We must therefore consider 

that they provide an approximate description. 

5.3. Continuous description 

Equation fl23l) is easily studied on a discrete lattice. However, when the wall position 
fiuctuates on the scale of many lattice distances (which it will do close to the critical line 
a = /3), we may replace the index /c by a continuous variable < A; < L. We will write 
Q{k,t) instead of Qkif) to indicate that we have performed this operation. Equation 
(1231) then becomes the Fokker-Planck equation 

dQ{k,t) dQ d'Q 

in which 

6 = D+-D^, D = \{D+ + D_) (26) 

appear as the drift velocity and the diffusion constant, respectively, of the domain wall. 
The boundary conditions become dQ{k,t)/dk = Q{k,t)/^ for k = 0,L. Subject to 
these, and if we set C, = D/5, the stationary solution of 1^5^ is 

k/£ 

^st(fc) = ^^J^,_^^ , 0<k<L. (27) 

For ^ positive (negative) the wall is localized near k = L (near /c = 0), as anticipated; 
moreover, the domain wall theory gives for the localization length |^| the explicit 
expression 

% That is, starting from the master equation that defines the time evolution of the particle 
configurations. 
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where the second expression was derived with the aid of (ITS]) . 

Let pL{k) denote the particle density at k averaged over the domain wall 
fluctuations; we will refer to this function as the 'density profile'. The domain wall 
theory gives for the density profile the following analytic expression, 

PL{k) =p' f dk' Q,,{k') + p'' [ dk'Q,,{k'). (29) 

Jo Jk 

Working this out with the aid of fl271) we find 

p,{k) =p' + (p^ - p') ^ _ ^„^^^ (30) 

which is valid for L ^ 1 and k of order L (as we used fl27j) ) and represents an exponential 
profile connecting p(0) = p^ at A; = to p(L) = p'^ ai k = L. We recall that 

P" = T^ , (31) 

1 + a 

/ = /, \ ,x (32) 



Equation (1281) shows that if we vary a (hence a) at fixed /3, the localization length |^| 
diverges at the critical value a = Oc = P- Exactly on the critical line a = P, the wall 
position distribution f l27|) becomes fiat, i.e., Qst{k) = l/L, the densities on each side 
of the wall are related by p'^ = p^/oc, and concomitantly the profile ( 15^ becomes the 
straight line 



PLik) 



Oc k 



ar 



(33) 



l + Oc 

where 1 — e""'' = ac- 

In figure |6] we show, for a system of length L = 302, for /3 fixed and for various 
values of a around Oc = /?, the theoretical density profiles piik) based on equations (128|) . 
( 130|) . (13T]) . and (!32|) . Also shown are the simulation data for piij^)- The error bars are 
of the order of the symbol sizes. For a = Uc the simulation results are indistinguishable 
from the theoretical straight line fl33|) . However, for a ^ a^ there is a small but definite 
discrepancy between the theory and the simulation data. We attribute this discrepancy 
to the approximate nature, mentioned above, of the master equation (!23l) on which the 
theory is based. Finding a theoretical description more accurate than (!23l) is at this 
moment an open challenge. 

5.4. Critical region 

Near criticality we may expand fl28l) in powers of Aa = a — Uc and obtain to leading 
order 

^ " JFAa ar Aa ^ ^ 
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Figure 6. Spatial density profiles as a function of the lattice site index A; in a 
system of length L = 300. The exit probability is /? = 0.4 and the entrance 
probabilities a are close to the critical value Oc = 0.4. The data points join 
together to produce heavy solid curves. A few error bars are given for the case 
a = 0.41. For all curves the error bars are of the order of the symbol sizes. 
Each simulation had a duration of 8 x 10^ time steps. Dashed curves: the 
theoretical predictions according to equations (f28]) . ([30l) . (f3T]) . and ([32]) . 



To lowest order in Aa we have p^ ~ p^/^c and thus the density profile fl30l) becomes 



PL{k) 



Otc. 



l + (ac 



-(L-k)/i 



-L/i 



^l + ac)(l-e- 



-{L-k)/i' 



-L/i 



(35) 



5.5. Scaling limit 

We shall now consider (135!) in the scaling limit Aa — )■ and L — ;■ oo with AaL fixed. 
In this limit the density profile piik) becomes a function of the single ratio k/L, as one 
may see by combining fl35|) and flMl) . and noticing that the variable z defined by 

ar Aa L 



lim — 

L— >-oo ^ 



1 + ttc 



at: 



(36) 
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Figure 7. Spatial density profiles as a function of the scaled lattice site index 
k/L, obtained from Monte Carlo simulations for various system sizes L and 
various a values, with {a—ac)L = —1.5 fixed. In the scaling limit the simulation 
data converge to a profile that is seen to depend only on k/L, as predicted by 
the domain wall theory. Dashed curve: the domain wall theoretical prediction 
for the density profile in the scaling limit. The exit probability /3 = 0.4 is fixed. 



has a fixed value. In this limit the profile (1351) may be simplified to 

1 - a, e^^^^^'^/-^) -e-'' 



Piik) 



1 + 



etc. 



(37) 



which depends only on the two scaling variables z and k/L. This scaling behavior of the 
profile is confirmed by the simulation data shown in figure [3 for increasing L the profiles 
converge to a limit curve. The limiting curve of the simulation data again exhibits a 
slight deviation from the theoretical curve ( 1371) . 

We return now to figure [5] and explain, again in the same scaling limit, the finite 
size effects observed there. In order to obtain the space averaged density we have to 
determine the mean value of the profile (137|) over < k/L < 1. The result is that 
we obtain for the L and Aa dependent space averaged particle density p the scaling 
expression 



P 



pi 



(Pf-Pc')^(^) 



(3^ 



Frozen shuffle update for open TASEP 19 

with — oo < 2; < oo and the scahng function \1/ given b}0 

Monte Carlo simulations were carried out on systems of different sizes L. We have 
seen in figure M that exactly at the critical point a = Oc the simulated profile is linear 
for all system sizes; as a consequence, in figure the space averaged density coincides 
for AaL = within error bars with the theoretical value. Within distance L~^ from 
criticality the transient time necessary to reach stationarity grows with L and we had to 
make sure that the data represent true stationary state values. The simulated system 
sizes L in figure E] are therefore smaller than the largest values shown in figure O The 
Monte Carlo data in figure [8] confirm the validity of the scaling as a function of the 
product variable {a — ac)L. 

As in figure [71 the Monte Carlo curves in figure M converge, when L increases, to 
a limit curve. Again, there is a slight but significant difference between the theoretical 
limit function and the simulation data, a consequence pointed out above of the fact that 
the domain wall theory is only approximate here. 

6. Conclusion 

Motivated by the modeling of pedestrian motion, we have in this paper pursued our 
study of a new 'frozen shuffle' update scheme for the totally asymmetric exclusion 
process (TASEP) on a one-dimensional lattice. An advantage of this scheme is that 
it fixes deterministic rules of priority when two hard-core particles attempt to hop to 
the same target site, a feature that allows for more efficient Monte-Carlo simulations. 
This is in contrast with the widely used parallel update which, when two pedestrians 
want to hop to the same site, generates conflicts which have to be solved by extra ad-hoc 
rules. Actually these algorithmic conflicts may have a physical counterpart: In dense 
crowds conflicts between real interacting pedestrians are thought to be one of the sources 
for the clogging effect observed near exits in emergency evacuation [3]. However, it is 
still an open question to determine whether the effect of the conflicts due to parallel 
update is overestimated or not, compared to the effect of real-life conflicts between 
pedestrians. Shuffle update yields an alternative that solves conflicts in a smoother way 
and may be more appropriate for pedestrian modelling in particular when the density 
is not too high. Besides, this scheme is easily accessible to analytical study. 

Compared to other update schemes, the frozen shuffle update modifies the 
properties of the stationary state in the closed system studied in reference [9] as well 
as in the open system studied here. The open system, just like the closed one, appears 
capable of sustaining 'free flow' and 'jammed' states. More specifically, the jammed 

+ The similarity of this expression to (|15p is probably a coincidence. 



Frozen shuffle update for open TASEP 



20 



0.8 



0.6 



0.4 



A 


L= 60 


O 


L= 15 


D 


L= 150 


* 


L = 240 


A 


L = 300 


V 


L=600 


— 


theory 




Figure 8. Space averaged particle density for various system sizes L, as a 
function of the scaling variable (a — ac)L. Solid curve: the theoretical scaling 
function according to equations ([38|) . (f39|) . and (f36|) . which describes the finite 
size rounding of the transition seen in figure [5j The Monte Carlo data are for 
/3 = 0.4, hence Oc = 0.4. For increasing L they converge to a scaling function 
whose deviation from the theoretical prediction is very small but nevertheless 
statistically significant. 



state's behavior is fully driven by the platoon structure, where platoons are ensembles 
of neighbor particles that move as a whole at each time step. As the platoon structure 
depends on the injection procedure, the flux of the jammed phase for the frozen shuffle 
update (in contrast to other update schemes) depends not only on the exit probability 
(3 but also on the entrance probability a. 

In spite of this specificity, the transition line between the free fiow and the jammed 
phase is found to be the usual a = /3 line. As we considered only the deterministic 
version of the model (particles move forward with probability 1 when the target site is 
empty), there is no maximal current phase. It could be of interest to extend the study 
to the case of particles hopping with probability p < 1. 

We have accompanied our theoretical results by Monte Carlo simulations, and 
found excellent agreement for the phase diagram. 

Finally, in order to explain the finite system size effects observed in the Monte Carlo 
simulations and especially important near the transition, we have applied a domain 
wall approach originally developed in references [IHl [161 [13 [El [19]. Frozen shuffle 
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update with deterministic dynamics leads to the remarkable property that the position 
of the domain wall that separates two homogeneous regions can be determined at the 
microscopic scale of the lattice distance. The domain wall theory appears to account 
quite well for the finite size effects, but we point out that nevertheless a small but definite 
discrepancy persists between this theory and the simulation data. The discrepancy 
concerns in particular the wall localization length, near the border of the system, for 
reasons still to be understood. 

An update scheme related to ours and preceding it, called 'random shuffle' update, 
has been proposed [5l [7] for the improvement of pedestrians models. That scheme has 
no obvious interpretation in terms of the interaction between pedestrians and should 
rather be viewed as an algorithm devised just to avoid conflicts. By contrast, the 'frozen 
shuffle' variant presented in this paper has a physical interpretation in terms of phase 
differences in the stepping cycle, which makes it relevant for applications to pedestrians. 

In forthcoming work |20j we will study the case of two intersecting lanes of 
pedestrians (TASEPs) simulated with frozen shuffle update. 
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